Question 1

Use R package UScensus2010county to complete the following tasks: (20 pt.)

if(!require(UScensus2010)) {install.packages("UScensus2010")}
Loading required package: UScensus2010
Loading required package: maptools
Checking rgeos availability: TRUE
Loading required package: foreign


Package UScensus2010: US Census 2010 Suite of R Packages
Version 0.11 created on 2011-11-18.

Zack Almquist, University of California-Irvine
ne

For citation information, type citation("UScensus2010").
Type help(package=UScensus2010) to get started.
library(UScensus2010)
if(!require(UScensus2010county)) {install.county("osx")}
Loading required package: UScensus2010county

UScensus2010county: US Census 2010 County Level Shapefiles and Additional
Demographic Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package="UScensus2010county") to get started.

For citation information, type citation("UScensus2010county").
library(UScensus2010county)
library(tmap)
package ‘tmap’ was built under R version 3.4.3

Question 1(a)

Plot a map of New York counties using the plot function.

data(new_york.county10)
shpfile <- new_york.county10
df1 <- shpfile@data
df1
plot(shpfile)

Question 1(b)

Plot a map of New York counties using the qtm function.

install.packages("tmap")
Error in install.packages : Updating loaded packages
library(tmap)
qtm(shpfile)

Question 1(c)

How many counties in New York State?

df
function (x, df1, df2, ncp, log = FALSE) 
{
    if (missing(ncp)) 
        .Call(C_df, x, df1, df2, log)
    else .Call(C_dnf, x, df1, df2, ncp, log)
}
<bytecode: 0x1074b8040>
<environment: namespace:stats>
print("There are 62 counties in New York State")
[1] "There are 62 counties in New York State"

Question 1(d)

What’s the 3-digit fips code of Broome County?

print("36007")
[1] "36007"

Question 1(e)

Compute descriptive statistics of the population column (P0010001), including total, minimum, maximum, mean, median, and skewness.

library(moments)
nypop <- df$P0010001
Error in df$P0010001 : object of type 'closure' is not subsettable

Question 1(f)

Create a histogram and a boxplot of the population.

rr a <- hist(nypop)

Error in hist.default(nypop) : 'x' must be numeric

Question 2

Use R package UScensus2010tract to complete the following tasks: (20 pt.)

rr if(!require(UScensus2010tract)) install.tract()

Loading required package: UScensus2010tract
there is no package called ‘UScensus2010tract’trying URL 'http://lakshmi.calit2.uci.edu/census2000/R/src/contrib/UScensus2010tract_1.00.tar.gz'
Content type 'application/x-gzip' length 190473496 bytes (181.6 MB)
==================================================
downloaded 181.6 MB

* installing *source* package ‘UScensus2010tract’ ...
** R
** data
** inst
** preparing package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded

UScensus2010tract: US Census 2010 Tract Level Shapefiles and Additional Demographic
Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package=\UScensus2010tract\) to get started.

For citation information, type citation(\UScensus2010tract\).
* DONE (UScensus2010tract)

The downloaded source packages are in
    ‘/private/var/folders/l3/pkmg9p591775dgyt86k3s9ym0000gn/T/RtmpIedmm3/downloaded_packages’
NULL

rr library(UScensus2010tract)


UScensus2010tract: US Census 2010 Tract Level Shapefiles and Additional Demographic
Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package=\UScensus2010tract\) to get started.

For citation information, type citation(\UScensus2010tract\).

Question 2(a)

Plot a map of New York census tracts using the plot function.

rr data(new_york.tract10) shpfile <- new_york.tract10 df <- shpfile@data df r plot(shpfile)

Question 2(b)

Compute the total population based on census tracts.

df
function (x, df1, df2, ncp, log = FALSE) 
{
    if (missing(ncp)) 
        .Call(C_df, x, df1, df2, log)
    else .Call(C_dnf, x, df1, df2, ncp, log)
}
<bytecode: 0x10d7b5578>
<environment: namespace:stats>
newpop <- df$P0010001
Error in df$P0010001 : object of type 'closure' is not subsettable

Question 2(c)

Select all census tracts in Broome County and plot the map.

broome <- county(fips="007",name="Broome",state="NY",level="tract")
Error in county(fips = "007", name = "Broome", state = "NY", level = "tract") : 
  could not find function "county"

Question 2(d)

What’s the total population of Broome County?

rr bpop <- sum(broome$P0010001) bpop

[1] 200600

Question 2(e)

Create a histogram and a boxplot of population based on census tracts of Broome County.

rr hist(broome$P0010001)

rr boxplot(broome$P0010001)

Question 2(f)

Select the first five columns of the shapefile of Broome County census tract; add a new population ratio column (= census tract population / county population); save the new shapefile to the hard drive.

popra <- data.frame(broome[,1:5])
Error in data.frame(broome[, 1:5]) : object 'broome' not found

Question 3

Use R packages raster and ncdf4 to complete the following tasks: (20 pt.)

library(raster)
library(ncdf4)

Question 3(a)

Load the multi-band raster dataset “NDVI.nc” into R.

ndvi <- raster("NDVI.nc")

Question 3(b)

Get the basic information about the dataset, including the number of rows, columns, cells, and bands; spatial resolution, extent, bounding box, and projection.

nrow(ndvi)
[1] 360
ncol(ndvi)
[1] 720
ncell(ndvi)
[1] 259200
nbands(ndvi)
[1] 36
res(ndvi)
[1] 0.5 0.5
extent(ndvi)
class       : Extent 
xmin        : -180 
xmax        : 180 
ymin        : -90 
ymax        : 90 
bbox(ndvi)
    min max
s1 -180 180
s2  -90  90
projection(ndvi)
[1] "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"

Question 3(c)

Aggregate all bands to generate a mean NDVI raster and a maximum NDVI raster; save the two new raster datasets to the hard drive.

ndvi.rb <- brick("ndvi.nc")
meanndvi <- mean(ndvi.rb)
plot(meanndvi)

writeRaster(meanndvi,filename = "meanndvi.tif",overwrite=TRUE)
maxndvi <- max(ndvi.rb)
plot(maxndvi)

writeRaster(maxndvi,filename = "maxndvi.tif",overwrite=TRUE)

Question 3(d)

Plot the maps of monthly NDVI of the year 2001

ndvi2k <- ndvi.rb[[13:24]]
plot(ndvi2k)

Question 3(e)

Create histograms of monthly NDVI of the year 2001.

hist(ndvi2k)
size changed to n because it cannot be larger than n when replace is FALSE39% of the raster cells were used (of which 42% were NA). 58093 values used.

Question 3(f)

Plot the NDVI map of July 2000; click any location with data on the map and retrieve the NDVI time series for all years; plot the NDVI time series of the selected location.

plot(ndvi.rb,7) 
location <- ndvi.rb[150, 50]
location <- click(ndvi.rb, n=1, xy=TRUE)
location <- click(ndvi.rb, n=1, xy=FALSE)

LS0tCnRpdGxlOiAiR2VvZzUzMyBMYWIgMTAiCmF1dGhvcjogIkpvc2h1YSBHb256YWxleiIKb3V0cHV0OgogIGh0bWxfZG9jdW1lbnQ6CiAgICB0b2M6IHllcwogIGh0bWxfbm90ZWJvb2s6CiAgICB0b2M6IHllcwogICAgdG9jX2Zsb2F0OiB5ZXMKLS0tCgojIyBRdWVzdGlvbiAxClVzZSBSIHBhY2thZ2UgVVNjZW5zdXMyMDEwY291bnR5IHRvIGNvbXBsZXRlIHRoZSBmb2xsb3dpbmcgdGFza3M6ICAoMjAgcHQuKQpgYGB7cn0KaWYoIXJlcXVpcmUoVVNjZW5zdXMyMDEwKSkge2luc3RhbGwucGFja2FnZXMoIlVTY2Vuc3VzMjAxMCIpfQpsaWJyYXJ5KFVTY2Vuc3VzMjAxMCkKaWYoIXJlcXVpcmUoVVNjZW5zdXMyMDEwY291bnR5KSkge2luc3RhbGwuY291bnR5KCJvc3giKX0KbGlicmFyeShVU2NlbnN1czIwMTBjb3VudHkpCmxpYnJhcnkodG1hcCkKYGBgCiMjIyBRdWVzdGlvbiAxKGEpClBsb3QgYSBtYXAgb2YgTmV3IFlvcmsgY291bnRpZXMgdXNpbmcgdGhlIHBsb3QgZnVuY3Rpb24uCmBgYHtyfQpkYXRhKG5ld195b3JrLmNvdW50eTEwKQpzaHBmaWxlIDwtIG5ld195b3JrLmNvdW50eTEwCmRmMSA8LSBzaHBmaWxlQGRhdGEKZGYxCnBsb3Qoc2hwZmlsZSkKYGBgCgojIyMgUXVlc3Rpb24gMShiKQkKUGxvdCBhIG1hcCBvZiBOZXcgWW9yayBjb3VudGllcyB1c2luZyB0aGUgcXRtIGZ1bmN0aW9uLgpgYGB7cn0KaW5zdGFsbC5wYWNrYWdlcygidG1hcCIpCmxpYnJhcnkodG1hcCkKcXRtKHNocGZpbGUpCmBgYAoKCiMjIyBRdWVzdGlvbiAxKGMpCQpIb3cgbWFueSBjb3VudGllcyBpbiBOZXcgWW9yayBTdGF0ZT8KYGBge3J9CmRmCnByaW50KCJUaGVyZSBhcmUgNjIgY291bnRpZXMgaW4gTmV3IFlvcmsgU3RhdGUiKQpgYGAKCiMjIyBRdWVzdGlvbiAxKGQpCQpXaGF04oCZcyB0aGUgMy1kaWdpdCBmaXBzIGNvZGUgb2YgQnJvb21lIENvdW50eT8KYGBge3J9CnByaW50KCIzNjAwNyIpCmBgYAoKIyMjIFF1ZXN0aW9uIDEoZSkJCkNvbXB1dGUgZGVzY3JpcHRpdmUgc3RhdGlzdGljcyBvZiB0aGUgcG9wdWxhdGlvbiBjb2x1bW4gKFAwMDEwMDAxKSwgaW5jbHVkaW5nIHRvdGFsLCBtaW5pbXVtLCBtYXhpbXVtLCBtZWFuLCBtZWRpYW4sIGFuZCBza2V3bmVzcy4gCmBgYHtyfQpsaWJyYXJ5KG1vbWVudHMpCm55cG9wIDwtIGRmJFAwMDEwMDAxCnN1bShueXBvcCkKbWluKG55cG9wKQptYXgobnlwb3ApCm1lYW4obnlwb3ApCm1lZGlhbihueXBvcCkKc2tld25lc3Mobnlwb3ApCmBgYAoKIyMjIFF1ZXN0aW9uIDEoZikJCkNyZWF0ZSBhIGhpc3RvZ3JhbSBhbmQgYSBib3hwbG90IG9mIHRoZSBwb3B1bGF0aW9uLgpgYGB7cn0KYSA8LSBoaXN0KG55cG9wKQphCmIgPC0gYm94cGxvdChueXBvcCkKYgpgYGAKCgojIyBRdWVzdGlvbiAyClVzZSBSIHBhY2thZ2UgVVNjZW5zdXMyMDEwdHJhY3QgdG8gY29tcGxldGUgdGhlIGZvbGxvd2luZyB0YXNrczogICAgKDIwIHB0LikKYGBge3J9CmxpYnJhcnkoVVNjZW5zdXMyMDEwKQppZighcmVxdWlyZShVU2NlbnN1czIwMTB0cmFjdCkpIGluc3RhbGwudHJhY3QoIm9zeCIpCmxpYnJhcnkoVVNjZW5zdXMyMDEwdHJhY3QpCmBgYAoKIyMjIFF1ZXN0aW9uIDIoYSkJClBsb3QgYSBtYXAgb2YgTmV3IFlvcmsgY2Vuc3VzIHRyYWN0cyB1c2luZyB0aGUgcGxvdCBmdW5jdGlvbi4KYGBge3J9CmRhdGEobmV3X3lvcmsudHJhY3QxMCkKc2hwZmlsZSA8LSBuZXdfeW9yay50cmFjdDEwCgpkZiA8LSBzaHBmaWxlQGRhdGEKZGYKcGxvdChzaHBmaWxlKQpgYGAKCiMjIyBRdWVzdGlvbiAyKGIpCkNvbXB1dGUgdGhlIHRvdGFsIHBvcHVsYXRpb24gYmFzZWQgb24gY2Vuc3VzIHRyYWN0cy4KYGBge3J9CmRmCm5ld3BvcCA8LSBkZiRQMDAxMDAwMQoKc25ldyA8LSBzdW0obmV3cG9wKQpzbmV3CmBgYAoKIyMjIFF1ZXN0aW9uIDIoYykKU2VsZWN0IGFsbCBjZW5zdXMgdHJhY3RzIGluIEJyb29tZSBDb3VudHkgYW5kIHBsb3QgdGhlIG1hcC4gCmBgYHtyfQpicm9vbWUgPC0gY291bnR5KGZpcHM9IjAwNyIsbmFtZT0iQnJvb21lIixzdGF0ZT0iTlkiLGxldmVsPSJ0cmFjdCIpCnBsb3QoYnJvb21lKQpgYGAKCiMjIyBRdWVzdGlvbiAyKGQpCldoYXTigJlzIHRoZSB0b3RhbCBwb3B1bGF0aW9uIG9mIEJyb29tZSBDb3VudHk/CmBgYHtyfQpicG9wIDwtIHN1bShicm9vbWUkUDAwMTAwMDEpCmJwb3AKYGBgCgojIyMgUXVlc3Rpb24gMihlKQpDcmVhdGUgYSBoaXN0b2dyYW0gYW5kIGEgYm94cGxvdCBvZiBwb3B1bGF0aW9uIGJhc2VkIG9uIGNlbnN1cyB0cmFjdHMgb2YgQnJvb21lIENvdW50eS4KYGBge3J9Cmhpc3QoYnJvb21lJFAwMDEwMDAxKQoKYm94cGxvdChicm9vbWUkUDAwMTAwMDEpCmBgYAoKIyMjIFF1ZXN0aW9uIDIoZikKU2VsZWN0IHRoZSBmaXJzdCBmaXZlIGNvbHVtbnMgb2YgdGhlIHNoYXBlZmlsZSBvZiBCcm9vbWUgQ291bnR5IGNlbnN1cyB0cmFjdDsgYWRkIGEgbmV3IHBvcHVsYXRpb24gcmF0aW8gY29sdW1uICg9IGNlbnN1cyB0cmFjdCBwb3B1bGF0aW9uIC8gY291bnR5IHBvcHVsYXRpb24pOyBzYXZlIHRoZSBuZXcgc2hhcGVmaWxlIHRvIHRoZSBoYXJkIGRyaXZlLiAKYGBge3J9CnBvcHJhIDwtIGRhdGEuZnJhbWUoYnJvb21lWywxOjVdKQpwb3ByYSRyYXRpbyA8LSBicm9vbWUkUDAwMTAwMDEvMjAwNjAwCnBvcHJhCmBgYAoKCiMjIFF1ZXN0aW9uIDMKClVzZSBSIHBhY2thZ2VzIHJhc3RlciBhbmQgbmNkZjQgdG8gY29tcGxldGUgdGhlIGZvbGxvd2luZyB0YXNrczogICAgICgyMCBwdC4pCmBgYHtyfQpsaWJyYXJ5KHJhc3RlcikKbGlicmFyeShuY2RmNCkKYGBgCgojIyMgUXVlc3Rpb24gMyhhKQkJCkxvYWQgdGhlIG11bHRpLWJhbmQgcmFzdGVyIGRhdGFzZXQg4oCcTkRWSS5uY+KAnSBpbnRvIFIuCmBgYHtyfQpuZHZpIDwtIHJhc3RlcigiTkRWSS5uYyIpCmBgYAoKIyMjIFF1ZXN0aW9uIDMoYikJCQpHZXQgdGhlIGJhc2ljIGluZm9ybWF0aW9uIGFib3V0IHRoZSBkYXRhc2V0LCBpbmNsdWRpbmcgdGhlIG51bWJlciBvZiByb3dzLCBjb2x1bW5zLCBjZWxscywgYW5kIGJhbmRzOyBzcGF0aWFsIHJlc29sdXRpb24sIGV4dGVudCwgYm91bmRpbmcgYm94LCBhbmQgcHJvamVjdGlvbi4KYGBge3J9Cm5yb3cobmR2aSkKbmNvbChuZHZpKQpuY2VsbChuZHZpKQpuYmFuZHMobmR2aSkKcmVzKG5kdmkpCmV4dGVudChuZHZpKQpiYm94KG5kdmkpCnByb2plY3Rpb24obmR2aSkKYGBgCgojIyMgUXVlc3Rpb24gMyhjKQkKQWdncmVnYXRlIGFsbCBiYW5kcyB0byBnZW5lcmF0ZSBhIG1lYW4gTkRWSSByYXN0ZXIgYW5kIGEgbWF4aW11bSBORFZJIHJhc3Rlcjsgc2F2ZSB0aGUgdHdvIG5ldyByYXN0ZXIgZGF0YXNldHMgdG8gdGhlIGhhcmQgZHJpdmUuIApgYGB7cn0KbmR2aS5yYiA8LSBicmljaygibmR2aS5uYyIpCm1lYW5uZHZpIDwtIG1lYW4obmR2aS5yYikKcGxvdChtZWFubmR2aSkKd3JpdGVSYXN0ZXIobWVhbm5kdmksZmlsZW5hbWUgPSAibWVhbm5kdmkudGlmIixvdmVyd3JpdGU9VFJVRSkKbWF4bmR2aSA8LSBtYXgobmR2aS5yYikKcGxvdChtYXhuZHZpKQp3cml0ZVJhc3RlcihtYXhuZHZpLGZpbGVuYW1lID0gIm1heG5kdmkudGlmIixvdmVyd3JpdGU9VFJVRSkKYGBgCgojIyMgUXVlc3Rpb24gMyhkKQkJCQpQbG90IHRoZSBtYXBzIG9mIG1vbnRobHkgTkRWSSBvZiB0aGUgeWVhciAyMDAxCmBgYHtyfQpuZHZpMmsgPC0gbmR2aS5yYltbMTM6MjRdXQpwbG90KG5kdmkyaykKYGBgCgojIyMgUXVlc3Rpb24gMyhlKQkKQ3JlYXRlIGhpc3RvZ3JhbXMgb2YgbW9udGhseSBORFZJIG9mIHRoZSB5ZWFyIDIwMDEuCmBgYHtyfQpoaXN0KG5kdmkyaykKYGBgCgojIyMgUXVlc3Rpb24gMyhmKQkJCQpQbG90IHRoZSBORFZJIG1hcCBvZiBKdWx5IDIwMDA7IGNsaWNrIGFueSBsb2NhdGlvbiB3aXRoIGRhdGEgb24gdGhlIG1hcCBhbmQgcmV0cmlldmUgdGhlIE5EVkkgdGltZSBzZXJpZXMgZm9yIGFsbCB5ZWFyczsgcGxvdCB0aGUgTkRWSSB0aW1lIHNlcmllcyBvZiB0aGUgc2VsZWN0ZWQgbG9jYXRpb24uIApgYGB7cn0KcGxvdChuZHZpLnJiLDcpIApsb2NhdGlvbiA8LSBuZHZpLnJiWzE1MCwgNTBdCmxvY2F0aW9uIDwtIGNsaWNrKG5kdmkucmIsIG49MSwgeHk9VFJVRSkKbG9jYXRpb24gPC0gY2xpY2sobmR2aS5yYiwgbj0xLCB4eT1GQUxTRSkKYGBgCgo=